#delimit ;
version 8.0;
*Harm Rival Interventions - Predicted Probabilities and Marginal Effects*
use "C:\Users\fosterdm\Downloads\KellerGrantFoster.dta";
nbreg _harm dempres unified election misery qtrgrw cap mid numstate monthsin LOC_ovr, robust;
preserve;
*Commanding run of Monte Carlo simulations to generate predicted probabilities and marginal effects - see Brambor, Clark, and Golder, 2005*
drawnorm b1-b12, n(10000) means(e(b)) cov(e(V)) clear;
postfile mypost predv_1 lo1 hi1 predv_2 lo2 hi2 predv_3 lo3 hi3 predv_4 lo4 hi4 predv_5 lo5 hi5 diff1_2 difflo1_2 diffhi1_2 diff2_3 difflo2_3 diffhi2_3 diff3_4 difflo3_4 diffhi3_4 diff4_5 difflo4_5 diffhi4_5 
using "harm3.tab", replace;
noisily display "start";
            
scalar h_dem=0;
scalar h_uni=0;
scalar h_elec=0;
scalar h_mis=9.585062;
scalar h_grw=.0174967;
scalar h_cap=.1896773;
scalar h_mid=12.1756;
scalar h_num=135.9286;
scalar h_mon=39.22321;
scalar h_loc=.3655;
scalar h_constant=1;

*Generating specific IV values of interest while holding values of controls constant (scalars above) - shifting presidential LOC from population min to -1SD to mean to +1SD to population max*
generate x_betahat1 = b1*h_dem+b2*h_uni+b3*h_elec+b4*h_mis+b5*h_grw+b6*h_cap+b7*h_mid+b8*h_num+b9*h_mon+b10*h_loc+b11*h_constant+b12*1;
generate x_betahat2 = b1*h_dem+b2*h_uni+b3*h_elec+b4*h_mis+b5*h_grw+b6*h_cap+b7*h_mid+b8*h_num+b9*h_mon+b10*.4128901+b11*h_constant+b12*1;
generate x_betahat3 = b1*h_dem+b2*h_uni+b3*h_elec+b4*h_mis+b5*h_grw+b6*h_cap+b7*h_mid+b8*h_num+b9*h_mon+b10*.4380731+b11*h_constant+b12*1;
generate x_betahat4 = b1*h_dem+b2*h_uni+b3*h_elec+b4*h_mis+b5*h_grw+b6*h_cap+b7*h_mid+b8*h_num+b9*h_mon+b10*.4632561+b11*h_constant+b12*1;
generate x_betahat5 = b1*h_dem+b2*h_uni+b3*h_elec+b4*h_mis+b5*h_grw+b6*h_cap+b7*h_mid+b8*h_num+b9*h_mon+b10*.4849886+b11*h_constant+b12*1;

*Specifying negative binomial transformation for predicted probabilities*
gen pr1=exp(x_betahat1);
gen pr2=exp(x_betahat2);
gen pr3=exp(x_betahat3);
gen pr4=exp(x_betahat4);
gen pr5=exp(x_betahat5);

*Specifying first differences for marginal effects*
gen d1_2=pr2-pr1;
gen d2_3=pr3-pr2;
gen d3_4=pr4-pr3;
gen d4_5=pr5-pr4;

*Specifying means for Monte Carlo*
egen p_v_1=mean(pr1);
egen p_v_2=mean(pr2);
egen p_v_3=mean(pr3);
egen p_v_4=mean(pr4);
egen p_v_5=mean(pr5);

egen d_f1_2=mean(d1_2);
egen d_f2_3=mean(d2_3);
egen d_f3_4=mean(d3_4);
egen d_f4_5=mean(d4_5);


tempname postfile mypost predv_1 lo1 hi1 predv_2 lo2 hi2 predv_3 lo3 hi3 predv_4 lo4 hi4 predv_5 lo5 hi5 diff1_2 difflo1_2 diffhi1_2 diff2_3 difflo2_3 diffhi2_3 diff3_4 difflo3_4 diffhi3_4 diff4_5 difflo4_5 diffhi4_5; 

*Generating confidence intervals for pp/me*
_pctile pr1, p(10,90);
scalar `lo1' = r(r1);
scalar `hi1' = r(r2);  
_pctile pr2, p(10,90);
scalar `lo2' = r(r1);
scalar `hi2' = r(r2);
_pctile pr3, p(10,90);
scalar `lo3' = r(r1);
scalar `hi3' = r(r2);  
_pctile pr4, p(10,90);
scalar `lo4' = r(r1);
scalar `hi4' = r(r2);
_pctile pr5, p(10,90);
scalar `lo5' = r(r1);
scalar `hi5' = r(r2);  
 

_pctile d1_2, p(10,90);
scalar `difflo1_2' = r(r1);
scalar `diffhi1_2' = r(r2);  
_pctile d2_3, p(10,90);
scalar `difflo2_3' = r(r1);
scalar `diffhi2_3' = r(r2);  
_pctile d3_4, p(10,90);
scalar `difflo3_4' = r(r1);
scalar `diffhi3_4' = r(r2);  
_pctile d4_5, p(10,90);
scalar `difflo4_5' = r(r1);
scalar `diffhi4_5' = r(r2);  

scalar `predv_1' = p_v_1;
scalar `predv_2' = p_v_2;
scalar `predv_3' = p_v_3;
scalar `predv_4' = p_v_4;
scalar `predv_5' = p_v_5;

scalar `diff1_2' = d_f1_2;
scalar `diff2_3' = d_f2_3;
scalar `diff3_4' = d_f3_4;
scalar `diff4_5' = d_f4_5;

post mypost (`predv_1') (`lo1') (`hi1') (`predv_2') (`lo2') (`hi2') (`predv_3') (`lo3') (`hi3') (`predv_4') (`lo4') (`hi4') (`predv_5') (`lo5') (`hi5')
 (`diff1_2') (`difflo1_2') (`diffhi1_2') (`diff2_3') (`difflo2_3') (`diffhi2_3') (`diff3_4') (`difflo3_4') (`diffhi3_4') (`diff4_5') (`difflo4_5') (`diffhi4_5');
     
*Estimates sent to post file; values manually transferred from post file 'harm3.tab' to Excel for graphing (see related Excel file)*
